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ABSTRACT 

Flux limited X-ray surveys of galaxy clusters show that clusters come in two roughly equally propor- 
tioned varieties: "cool core" clusters (CCs) and non-"cool core" clusters (NCCs). In previous work, 
we have demonstrated using cosmological A-body + Eulerian hydrodynamic simulations that NCCs 
are often consistent with early major mergers events that destroy embryonic CCs. In this paper we 
extend those results and conduct a series of simulations using different methods of gas cooling, and 
of energy and metal feedback from supernovae, where we attempt to produce a population of clusters 
with realistic central cooling times, entropies, and temperatures. We find that the use of metallicity- 
dependent gas cooling is essential to prevent early overcooling, and that adjusting the amount of 
energy and metal feedback can have a significant impact on observable X-ray quantities of the gas. 
We are able to produce clusters with more realistic central observable quantities than have previously 
been attained. However, there are still significant discrepancies between the simulated clusters and 
observations, which indicates that a different approach to simulating galaxies in clusters is needed. 
We conclude by looking towards a promising subgrid method of modeling galaxy feedback in clusters 
which may help to ameliorate the discrepancies between simulations and observations. 
Subject headings: cosmology: theory — galaxies: clusters: general — hydrodynamics intergalactic 
medium — methods: numerical 



1. INTRODUCTION 

A wide variety of X-ray surveys (e.g., Sanderson et al. 
2006; Chen et al. 2007; Johnson ct al. 2009) of galaxy 
clusters give the interesting result that the proportion 
of "cool core" (CC) clusters and non-"cool core" (NCC) 
clusters is roughly equal. While authors use different 
specific definitions for CC clusters, they are fundamen- 
tally differentiated from NCC clusters by bright, cuspy 
central X-ray profiles indicative of high central densities 
and correspondingly short cooling times (< Hq 1 ), and 
low central temperatures (« 30%-40% of the virial tem- 
perature; Ikebe et al. 1997; Lewis et al. 2002; Peter- 
son et al. 2003). Early observations of CC clusters (Lea 
et al. 1973; Cowie & Binney 1977; Fabian & Nulsen 1977; 
Mathews & Bregman 1978) motivated the creation of the 
"cooling flow" model (see Fabian (1994) for a review) 
to explain the existence of CCs within galaxy clusters, 
wherein central gas quasi-hydrostatically cooled to low 
temperatures, and was replenished with hot gas flowing 
in from the intracluster medium (ICM). However, the 
cooling flow model suggests a higher star formation rate 
than is observed (McNamara & O'Connell 1989; Cardiol 
et al. 1998; Edge 2001), predicts a high mass deposition 
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rate (Makishima et al. 2001), and cooler central temper- 
atures (Peterson et al. 2001; Sanders et al. 2008), and 
therefore newer ideas have been created to explain the 
bimodality. 

Theories to explain the discrepancies often include ad- 
ditional sources and methods of transport of heat in a 
cluster. These ideas include thermal heat conduction 
(Zakamska & Narayan 2003; Voit 2011) combined with 
sound waves (Ruszkowski et al. 2004) or with turbulence 
(Dennis & Chandran 2005), conduction combined with 
cosmic rays (Loewenstein et al. 1991; Guo & Oh 2008), 
and active galactic nuclei (AGN) (Binney & Tabor 1995; 
Rephaeli & Silk 1995; McNamara & Nulsen 2007). Nu- 
merical simulations have demonstrated that AGN may 
be important sources of feedback energy and significantly 
impact the energy balance in clusters. Simulations show 
that AGN bring simulations closer to observations by 
preventing overcooling, which lowers cluster central den- 
sities (Sijacki et al. 2007; Teyssier et al. 2011), improves 
the distribution of metals (Fabjan ct al. 2010; McCarthy 
ct al. 2010), and regulates star formation rates and trig- 
gers quenching (Di Matteo et al. 2005; Sijacki & Springcl 
2006). 

Other phenomena that may explain the observed clus- 
ter population properties include those methods that dis- 
turb the gas at the centers of clusters, or which disrupt 
the cool core. This class of phenomena is most natu- 
rally explained by merger events between halos. It is, of 
course, impossible to observationally follow cluster merg- 
ers from start to finish and measure its effect on cluster 
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cores, but there is some evidence that NCC clusters are 
associated with past merger events. By measuring the 
offset between the brightest cluster galaxy and the X-ray 
centroid, Sanderson et al. (2009) find that clusters with 
larger offsets, which imply a state of dynamical distur- 
bance, have weaker CCs. In their sample, Rossetti et al. 
(2011) find no clusters with giant radio halos (which are 
associated with past merger events) that can be classi- 
fied as CC. Using detailed density profiles of observed 
CC and NCC clusters, Eckert ct al. (2012) find that as 
compared to CC, the outskirts of NCC clusters have flat- 
ter profiles, which they argue is evidence of past major 
merger events which redistribute gas between all regions 
of the cluster more effectively than centrally-dominated 
sources of entropy, such as AGN. Hcnning et al. (2009) 
find that observed NCC clusters are warmer in the pe- 
riphery than CC clusters, and that an analogous popula- 
tion of NCC clusters sampled from numerical simulations 
have a richer history of early merger events compared to 
their simulated CC cluster counterparts. 

Previously, our group has used cosmological simula- 
tions to explore the idea that mergers influence the cre- 
ation of NCC clusters. In Burns et al. (2008) (hereafter 
B08) we performed simulations that included gas cooling 
plus star formation with supernovae feedback, and we 
found that galaxy clusters are born cool, but may warm 
up when impacted by a major merger at some point in 
their early histories. Interestingly, Rossetti & Molcndi 
(2010) find most of their sampled observed NCC clus- 
ters have regions of high mctallicity and low entropy gas, 
which they interpret as indicators of a past cool phase for 
NCC clusters. Using simulations, ZuHone et al. (2010) 
show that when the core of a cluster interacts with a 
merging subcluster, the core gas can enter a "sloshing" 
phase that can delay the development of the cool core for 
one or more Gyr. In a related paper, ZuHone (2011) finds 
that cluster mergers over a wide range of mass ratios and 
impact parameters act to mix up the phases of gas and 
raise the overall entropy floor, making CC formation (or 
preservation) more difficult. 

The assumption that the galaxy cluster population is 
approximately evenly split between cool core and non- 
cool core clusters is likely too simplistic. It is proba- 
ble that there is an observational bias in the CC/NCC 
cluster ratio in surveys because flux-limited X-ray obser- 
vations will naturally over-sample the centrally brighter 
CC groups. Observations by Eckert et al. (2011) sug- 
gest that there is a bias towards CC clusters of as much 
as 30%. Hudson et al. (2010) (hereafter H10) used HI- 
FLUGCS (Highest X-ray FLUx Galaxy Clusters Sam- 
ple) data (Rciprich & Bohringcr 2002) to do a statistical 
analysis on a large number of cluster measureables in an 
attempt to establish outlines between the CC and NCC 
subgroups. For a few of the observables, most notably 
central cooling time and entropy (see their Fig. 4), they 
find that the best statistical fit in fact divides the popu- 
lation into three groups: NCCs, Weak CCs (WCCs), and 
Strong CCs (SCCs). We discuss the work presented in 
H10 in more detail in §3. 

In B08 we defined a CC cluster as one with a > 20% 
drop in central gas temperature compared to the sur- 
rounding gas. Using this definition, we produced both 
CC and NCC clusters in a single cosmological simula- 
tion for the first time. The work presented in this paper 



builds on the results of B08 in two ways. First, we focus 
on comparing our results with observations to understand 
the effects of the different physical models, and use these 
comparisons to modify our simulations in an effort to 
make more realistic cluster populations. Second, we in- 
clude two physically-motivated numerical models, one for 
gas cooling and one for metal and energy feedback, that 
we use in an attempt to address shortcomings in our pre- 
vious work. We compare and contrast our simulations 
in order to understand how the physical models affect 
the clusters, in a way inspired by the pioneering work of 
Lewis et al. (2000) that examined in detail the effect of 
gas cooling on simulated galaxy clusters. Many groups 
have investigated the simulation of clusters with gas cool- 
ing and supernovae feedback (e.g. Valdarnini 2003; Bor- 
gani 2004; Kravtsov et al. 2005; Kay et al. 2007; Nagai 
et al. 2007; Tornatore et al. 2007; Fabjan et al. 2011), but 
in this paper, we use similar physics - however we use a 
large set of clusters (~ 65 — 70 per simulation) and focus 
on only the properties of core of the clusters. As we will 
show, when additional cluster observables are considered, 
the B08 simulations fail to reproduce other observed clus- 
ter characteristics such as the central entropy. Therefore, 
our motivations for including the new models of cooling 
and feedback are to address clear deficits in the methods 
used in B08, improve the concordance with observations 
for our simulated clusters, and discover the shortcomings 
of current methods and how they might motivate future 
numerical models. 

We note that the simulations discussed in this paper 
do not include a prescription for AGN formation and 
feedback. We have already listed some of the benefits 
of AGN in cluster simulations, and it is well known that 
simulations that include only cooling + star formation 
do not do a good job of reproducing all of the observed 
properties of clusters (Sijacki & Springel 2006; Sijacki 
et al. 2007; Puchwein et al. 2008). This paper is, there- 
fore, an exploration of the limitations of models that use 
only stellar feedback, and quantifies how well (or poorly) 
this type of simulation performs when compared against 
standard cluster observable quantities. 

In Section 2 we describe our simulations and numer- 
ical methods, including the details of the gas cooling 
and feedback mechanisms. We base our comparisons to 
observations using the results presented in H10 of HI- 
FLUGCS data, which we describe in detail in §3. Our 
first four simulations, presented in §4, investigate the ef- 
fects of a variety of methods for gas cooling and stellar 
feedback between simulations, and we also compare to 
observations. In response to the results of §4, we per- 
form three additional simulations and compare them to 
observations in §5. In §6 we summarize the results and 
shortcomings of our simulations, and in §7 we present 
our conclusions and a discussion regarding the steps that 
must be taken in order to produce more realistic clusters. 

2. NUMERICAL METHODS 

Our simulations are performed using Enzo 4 (O'Shea 
et al. 2005a; Norman et al. 2007), an open-source, 
community-developed code that solves hydrodynamics 

4 http://enzo-project.org/; our simulations were run using 
the Enzo development branch with Mercurial revision identifier 
dl2ea971621c. 
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on an Eulerian mesh. In these simulations, we use the 
ZEUS scheme for solving hydrodynamics (Stone Sz Nor- 
man 1992) , with an TV-body scheme for evolving the dark 
matter and star particles. Enzo includes adaptive mesh 
refinement (AMR) capability, which wc discuss in more 
detail below. Enzo has been compared to other cos- 
mological codes with favorable results for ICM physics 
(Frenk et al. 1999; O'Shea et al. 2005b; Agertz et a!. 2007; 
Tasker et al. 2008; Vazza et al. 2011). We use a ACDM 
cosmology, and our initial conditions are generated at 
z = 99 using the CDM transfer function from Eiscn- 
stein & Hu (1999) using the parameters f2 m = 0.268, 
n A = 0.732, Q b = 0.0441, h= 0.704, cr 8 = 0.82, and 
n s = 0.97 (WMAP3; Spergel et al. 2007). 

Our simulations employ a 128 3 (h^ 1 Mpc) 3 volume 
with 256 3 top-level grid zones and dark matter parti- 
cles. Dark matter particles are given a mass of 7.8 x 10 9 
h _1 MQ, and the mean baryon mass of a root-grid cell is 
3.1 x 10 9 h~ 1 Mq . AMR is enabled throughout the entire 
volume on cells containing baryonic and/or dark matter 
density 8.0 times the mean at that level, for up to 5 ad- 
ditional levels. This results in a peak co-moving spatial 
resolution of 15.6 ft -1 kpc, which is sufficient to resolve 
features having a scale of ~ 50— 100 h^ 1 kpc such as the 
cool cores of clusters, but is insufficient to resolve any 
more detailed structure. 

2.1. Cooling Physics 

Radiative cooling of gas is applied on every cell at ev- 
ery time step. Using one of the two cooling methods de- 
scribed below, an amount of energy to be radiated away 
is calculated, and that energy is subtracted from the cell. 
For all simulations using the metal-dependent Cloudy 
cooling method (described below) , a uniform, metagalac- 
tic ionizing UV background (q a = 1.5, Haardt & Madau 
1996) from quasars is applied after z=7. 

Two of our simulations use metallicity-independent 
Raymond-Smith cooling (hereafter "RS cooling"; Ray- 
mond et al. 1976), the identical method to the one used in 
the B08 simulations. RS cooling rates are computed from 
an optically thin fully-ionized plasma emission model 
(Brickhouse et al. 1995; Sarazin & White 1987) that as- 
sumes a constant metallicity of 0.5 relative to solar (Zq). 
Analytical approximations of cooling rates are stored in a 
lookup table as a function of temperature, and the cool- 
ing rate of a cell is found by multiplying by the density 
squared of the cell. We have run a simulation with RS 
cooling using rate tables computed assuming a constant 
Z/Zq = 0.25 metallicity and did not find any significant 
changes from Z/Zq = 0.5. Using a cooling function de- 
pendent only on temperature is simplistic - in particular, 
the assumption of constant metallicity is questionable. 
As we will show in Section 4.1, the assumption of a con- 
stant, high metallicity leads to overcooling of gas at high 
redshift and substantially affects cluster observables. 

To address this issue, we replace RS cooling with a 
cooling model partially based on the photoionization 
code Cloudy (Ferland et al. 1998). Throughout this pa- 
per, we use the term "Cloudy cooling" to describe this 
model that includes metallicity- dependent cooling rates 
(Smith et al. 2008, 2011). The simulations that utilize 
Cloudy cooling actually combine two methods to calcu- 
late cooling rates, and the total rate is the sum of the 



two. For atomic H and He, the individual abundances are 
tracked and the non-equilibrium cooling rates are com- 
puted directly via a network of coupled equations (Abel 
ct al. 1997; Anninos et al. 1997). The metal cooling rates 
from all atomic species between Li and Zn (assuming so- 
lar abundances for relative species fractions), and many 
molecular species as well, are found by referencing pre- 
computed tables that depend on temperature, density, 
electron fraction, and (crucially) metallicity up to z=7, 
and also the UV background after z=7. The tables are 
built by using the Cloudy code to compute the equilib- 
rium rates of cooling over a grid of values in the phase 
space of the aforementioned dependencies, and the rates 
are interpolated from the grid of values by Enzo as the 
simulation runs. The range of the dependencies used to 
build the tables covers all the relevant physical values for 
our simulations. 

2.2. Star Formation and Feedback 

Just as radiative cooling is an important part of the 
energetics of a galaxy cluster, so is star formation and 
the associated feedback. As in B08, we include a pre- 
scription for star formation in our simulation using the 
formulation described in Cen & Ostriker (1992), which 
we briefly describe here. In all grid cells that are locally 
at the most highly-refined level, a collisionless "star par- 
ticle" is created in a cell if the following conditions are 
met in that cell: 1) the gas density is at least 100 times 
higher than the mean gas density in the volume; 2) the 
flow of gas is locally converging; 3) the gas cooling time is 
less than the total (gas + dark matter) dynamical time; 
and, 4) and the mass of gas in the cell exceeds the Jean's 
mass. The mass of the star M sta r = -^ceiiAi/idyn de- 
pends on the mass of gas in the cell M ce \\, the length of 
the time step At, and the dynamical time of the gas in 
the cell tdyn, which is calculated for every cell and has 
a minimum of 1 Myr. A star particle is formed only if 
the minimum mass constraint is satisfied (due to compu- 
tational limitations), which we set to M star > 10 9 M Q . 
The star particle is given a metallicity Z star where the 
value is same as the metallicity of the gas from which 
the particle is formed. 

After the star is formed, energy is deposited back into 
the gas to simulate feedback from Type II supernovae. 
The amount of energy feedback 

e = e E M star c 2 (1) 

is controlled by the dimensionless input parameter €e *C 
1. In all our simulations star particles also return met- 
als to the gas in order to model enrichment from super- 
novae. Recall that RS cooling assumes a constant metal- 
licity, and although star particles increase the metallic- 
ity of the surrounding gas in simulations employing RS 
cooling, the calculated cooling rates are completely inde- 
pendent of this effect. However, we will show that gas 
metallicity is very important in calculations that use the 
metallicity-dependent Cloudy cooling model. The total 
mass of metals Mz returned from the star particle to the 
gas phase is described by the equation 

M z = M star [(1 - Z stai )e z + 0.25Z star ] , (2) 

where ez < 1 is a dimensionless input to our simulations. 
Both the energy and metal feedback is applied over a 
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12idyn time period after the particle is formed, where 
the rates rise linearly for t < t^ yn and then decrease ex- 
ponentially after idyn- The total feedback quantities are 
independent of simulation time step, but the increased 
heat in the gas may shorten the time steps due to the 
Courant condition. The value of t^ yn varies for each star 
particle (it is assigned at creation), and in the cluster sim- 
ulations discussed in this paper it ranges from roughly 10 
to 30 Myr. See Table 1 for the values of ce and ez used 
in each of our simulations. 

In most grid-based simulations that include star forma- 
tion with feedback (including those in B08), energy and 
metals are deposited entirely within the cell that contains 
a particular star. As discussed in Smith et al. (2011), 
dumping all of the thermal or kinetic feedback into a 
single high density cell may result in calculated cooling 
times significantly shorter than a hydrodynamical time 
step, which results in overcooling the gas in the cell. In 
particular, overcooled gas prevents the stellar feedback 
from being spread out into the intracluster medium. In 
order to address this, Smith et al. implemented in Enzo 
the "distributed feedback" method, wherein energy and 
metal feedback are deposited into more than one cell sur- 
rounding the star particle. The same total quantity of 
feedback is injected into the gas regardless of the num- 
ber of cells, and each cell receives an equal share. Smith 
et al. show (e.g., their Fig. 2) that simulations that 
use Cloudy cooling + distributed feedback produce star 
formation rates that more closely resemble the mean cos- 
mic star formation history than those that use single-cell 
feedback. For all of our simulations that use it (see Ta- 
ble 1), distributed feedback is applied over a 3 x 3 x 3 
cube of cells at the highest resolution centered on the cell 
containing the star. 

2.3. Analysis Methods 

All simulations are evolved to and analyzed at z=0. 
Our analyses are performed using the yt 5 toolkit (Turk 
et al. 2011). Cluster dark matter halos are located using 
a parallel version of the HOP halo finder (Skory et al. 
2010), and only virialized clusters with total (gas + dark 
matter) mass M 2 oo greater than 10 14 M are kept for 
our samples. We define M200 as the mass enclosed inside 
a sphere centered on the cluster with average density 200 
times the mean background density. 

Central quantities and radial profiles of observed clus- 
ters are typically found first by centering on the peak of 
X-ray emission. X-ray emission is roughly proportional 
to p^T 1 / 2 , therefore we center our profiles on a gas den- 
sity peak close to the center of the cluster, but not neces- 
sarily the absolute gas density peak, for the following rea- 
son. In many of our simulated clusters, infalling clumps 
of gas are dense enough such that their cold centers are 
not shock heated immediately as they pass through the 
outer virialized gas of the cluster. This means that the 
absolute point of maximum gas density can be offset from 
the actual central X-ray peak of the cluster. The proper 
center of a cluster can be located by eye. However, this 
is tedious and slow, so we have developed an automated 
method that in our tests is nearly always (better than 
99%) identical to manual identification. 

5 http://yt-project.org/ 



First, we find gas density "clumps" inside the cluster 
defined by regions of gas of at least 0.025 times the maxi- 
mum gas density in the cluster. Within a factor of a few, 
the results are not very sensitive to the choice of density 
threshold, but values much higher or lower will result in 
incorrect identification of cluster centers. The clumps de- 
fined by the threshold are not necessarily gravitationally 
bound; they are simply defined by the volume enclosed 
by the density contour. Next, we identify the location 
and value of maximum density inside each clump. For 
smooth clusters, the point of absolute maximum density 
is coincident with the maximum inside the single clump. 
For clusters with significant substructure, however, there 
are multiple discrete clumps, and each has its own max- 
imum density point. To choose the proper X-ray peak, 
we pick the maximum density point of the clump that 
minimizes the equation 

Q—( ^ \ l V dustcr — V c i ump | 2 /MclustciA ^ 
\R — dJ |v c i ustc r| V M=lump / ' 

where R is the radius of the cluster, d is the distance 
from the center of dark matter mass of the cluster to 
the maximum density point of that clump, v c i us t cr is the 
bulk velocity of the cluster, v c i ump is the bulk velocity 
of the clump, M c i ump is the total mass of the clump (gas 
+ dark matter), and M c i ustOT is the total mass of the 
cluster. The first term enforces a preference for clumps 
closer to the overall center of mass, the second attempts 
to eliminate clumps of gas moving quickly compared to 
the bulk of the cluster (such as infalling clumps of gas), 
and the third chooses larger clumps over smaller (which 
has the effect of eliminating small clumps of infalling gas 
and is complementary to the second term). 

3. COMPARISON OBSERVATIONAL DATASET: HIFLUGCS 

Our main comparisons against observations use the 
set of clusters analyzed in H10, which includes all in 
the HIFLUGCS sample (originally defined in Rciprich & 
Bohringer 2002). The HIFLUGCS is a statistically com- 
plete, flux-limited sample of clusters with So.5-2/seV > 
2 x 10~ n ergs/sec/cm 2 . This sample is comprised of 64 
of the X-ray brightest clusters at high galactic latitudes. 
The average redshift for the rich clusters in the sample 
is 0.053. All clusters in the HIFLUGCS sample have 
been observed by both Chandra and XMM-Newton with 
good signal/noise, and most have been observed mul- 
tiple times. The median Chandra and XMM-Newton 
integration times are 65 ksec and 68 ksec, respectively. 
HIFLUGCS has a significant number of both cool core 
and non-cool core clusters, mostly Abell clusters, but 
also includes a sample of galaxy groups, as well. The HI- 
FLUGCS observations are analyzed in H10 by a CIAO 
+ CALDB pipeline described in detail in Hudson et al. 
(2006). The pipeline produces mosaiced X-ray images 
that are used to calculate 16 observable quantities for 
each cluster. 

The goal of work presented in H10 is to investigate 
which physical properties of clusters can be used, and 
with what confidence, to differentiate between types of 
clusters. Using a statistical test for Gaussian bimodality 
(or trimodality if the algorithm gave it a higher confi- 
dence; Ashman et al. 1994) the clusters are divided into 
two (or three) subgroups. 
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TABLE 1 

Simulation Parameters 



Simulation Label Cooling # of Feedback Cells t£ a ez h # of Clusters 

RS-Single RS 1 1.0 x 10~ 5 0.1 78 

RS-Dist RS 27 1.0 x lO" 5 0.1 79 

Cloudy-Single Cloudy 1 1.0 X lO" 5 0.1 67 

Cloudy-HZ-LE Cloudy 27 1.0 X 10~ 5 0.1 78 

Cloudy-LZ-LE Cloudy 27 1.0 X 10~ 5 0.02 79 

Cloudy-LZ-ME Cloudy 27 8.0 X 10~ 5 0.02 76 

Cloudy-LZ-HE Cloudy 27 2.0 X 10" 4 0.02 71 

Adiabatic N/A N/A N/A N/A 80 

a Dimensionless energy feedback parameter; e — e£;M s t ar c 2 . 

k Dimensionless metallicity feedback parameter; Mz — A/ Star ((l — Z st&T )ez + 0.25Z star ). 



The subgroups for each observable are labeled accord- 
ing to the physical interpretation of the differences be- 
tween subgroups. In the case of a bimodal distribution, 
the clusters are divided into the usual CC and NCC sub- 
groups. For trimodal distributions, the CC subgroup 
is divided into Strong-Cool Core (SCC) and Weak-Cool 
Core (WCC) subgroups. The cut(s) between subgroups 
are defined by the Gaussian test, and clusters are as- 
signed to subgroups based on where they fall relative to 
the cuts. 

The main results of H10 are summarized with his- 
tograms and Gaussian fits for each quantity in their Fig. 
4. All quantities show high confidence of at least bi- 
modality. Two quantities, central cooling time and en- 
tropy, are further deemed to be most likely trimodal. 
In order to establish the "defining" subgroups (those to 
which the Gaussian test gives the highest statistical like- 
lihood) between types of clusters, they compare the like- 
lihood of all bimodal distributions (including a bimodal 
distribution of central cooling time and entropy). They 
find that using bimodal distributions of cooling time and 
entropy give the most statistically likely subgroups, and 
therefore they describe these two as the defining param- 
eters of cluster type. As previously mentioned, they find 
that further dividing central cooling time and central en- 
tropy into three subgroups is actually more likely than 
two, resulting in the trimodal central cooling time and 
entropy subgroups as the defining subgroups for cluster 
type. The trimodal cuts for cooling time as defined by 
H10 are SCC < 1 h- 1 / 2 Gyr, 1 h,- 1 / 2 Gyr < WCC < 7.7 
h- 1 / 2 Gyr, and 7.7 h^ 1 ' 2 Gyr < NCC. For entropy they 
are SCC < 22 /i" 1 / 3 keV cm 2 , 22 hr 1 ! 3 keV cm 2 < WCC 
< 150 hr 1 ! 3 keV cm 2 , and 150 hr 1 ! 3 keV cm 2 < NCC. 
Most clusters lie in matching central cooling time and 
central entropy subgroups, but 5 clusters lie just across 
dividing lines and have different designations (see Fig. 5 
ofHIO). 

In B08 we define a CC cluster as one with a > 20% 
drop in temperature to the center when compared to the 
temperature of the gas where the drop begins (i.e. where 
the inward radial temperature profile slope becomes neg- 
ative). This is qualitatively similar to the definition used 
in H10, which we will adopt in this paper. They de- 
fine the central temperature ratio 6 as Xo/T v i r , where Tq 
is the central temperature, and T v ; r is the "temperature 

6 H10 uses the term "central temperature drop," which we mod- 
ify in recognition of the many observed clusters with To/T v ; r > 1. 



of the X-ray emitting gas that is in hydrostatic equilib- 
rium with the cluster potential." T v ir is found using a 
mass-scaling relation (see §3.2.4 for details). They find 
that a single bimodal cut between CC/NCC clusters of 
To/T v i T = 0.7 gives them one of their most statistically 
confident bimodal distributions. 

Due to the high confidence of these three parameters 
(central cooling time, entropy, and temperature ratio) 
as defining characteristics, and our previous use of one 
of them, we will focus our observational comparisons on 
these quantities. 

We make use of several other observational datasets for 
additional comparisons, but these are not used to con- 
strain our simulations. We compare our simulated cluster 
mctallicities to results presented in Matsushita (2011), 
and our star formation rates to values from Bouwens 
et al. (2007) and Hopkins (2004, 2007). 

3.1. The Definition of Central Region 

The definition of the central region in H10 is the vol- 
ume enclosed inside the innermost annular bin of the ra- 
dial profiles derived from the X-ray maps. The physical 
size of this region depends on many things, in particular 
the size of the cluster itself and the number of source 
counts. We cannot use this method for our simulated 
data - we can, however, infer an upper limit. The cen- 
tral cooling time is calculated by H10 (see their Eqn. 
15, and our §3.2.2) using the average temperature inside 
r < 0.048 i?5oo, which implies that the central region 
is no larger than this. For some of the smaller simu- 
lated clusters, using this radius for the central region 
would result in calculating quantities over just a few cells. 
To avoid stochasticity due to the inclusion of very small 
numbers of cells, we define the central region as the inner 
r < 50 kpc for all clusters, which covers ~ 50 maximum- 
resolution cells. This radius was chosen after investigat- 
ing the values of our observables over a range of central 
region radii. Using a radius less than 50 kpc gives a scat- 
ter of values that covers many orders of magnitude, and is 
much larger scatter than observations. Choosing a radius 
larger than 50 kpc showed very little change in cooling 
time beyond the 50 kpc values; no larger than a factor 
of two out to several hundred kpc. 50 kpc was therefore 
chosen as a compromise between volumes too small and 
volumes much larger than used in the observational data 
analyses. In fact, for the majority of the clusters, 50 kpc 
is within a factor of two of 0.048 -R500, and is therefore a 
reasonable estimate for the central radius. 
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3.2. Our Three Main Observables 

In the previous section, we identified the central cool- 
ing time, entropy, and temperature ratio as the three 
quantities upon which we will focus our comparisons. In 
this section, we describe how we calculate each observ- 
able quantity. 

3.2.1. Spectroscopic-Like Temperature 

In place of the simulated gas temperature T, we use 
the spectroscopic- like temperature T spoc , calculated as 
follows: 



J p 2 T ( a -0.5) d y 



(4) 



over the volume of interest, where a — 0.75. Previ- 
ous work has shown that spectroscopic-like temperature 
(see Mazzotta et al. 2004; Rasia et al. 2005) is a reason- 
able proxy for the measured X-ray spectral temperature. 
When making mock temperature maps, we perform a 
similar calculation, but use this weighting for the line- 
of-sight integral, as has been used in earlier studies (e.g., 
Hallman et al. 2010). This weighting has been shown 
to reproduce the fitted spectral temperature better than 
standard emission weighting or mass weighting of the 
temperature in simulations (Mazzotta et al. 2004). 

We have discovered in prior work that T spec has some 
critical limitations in the context of simulations with ra- 
diative cooling. T spec in its original incarnation is a cali- 
brated weighting, using simulations of clusters with mean 
temperatures greater than T — 2keV. Additionally, 
when calculating T spec , Mazzotta et al. (2004) ignored 
gas particles with T < O.bkeV. In regions in some of 
our simulated clusters where cooling dominates, the gas 
reaches much lower temperatures than this. Because gas 
at these low temperatures contributes negligibly to the 
X-ray emission, including these zones in the T spec calcu- 
lation results in temperatures that are not representative 
of the X-ray temperature for these clusters. Therefore, 
in order to more accurately model the X-ray temperature 
for the clusters, the cold gas should be removed from the 
calculation. We discuss the implications of this and our 
choice to include or not include the cold gas in clusters 
in §4.5. 

3.2.2. Central Cooling Time 

We calculate central cooling time to, cool over the central 
region Vq as follows (adapted from Eqn. 15 of H10): 



to 



cool 



3 K + n^f Vo kT spec dV 
2 C nf L A(T spec )dV 



(5) 



where n' e — \fjn^) is the number density of electrons (a 
volume- weighted average), v! i — n' e /( is the number den- 
sity of ions, £ is the ratio of electrons to ions and is set 
to 1.2 following H10, and A(T spoc ) is the cooling func- 
tion. Values of A come from the same cooling function 
table used in the RS simulations, which is very similar 
to the APEC (Smith et al. 2001) cooling function used 
in H10 for optically thin plasmas. Using a metallicity- 
independent cooling function for the cooling time cal- 
culation is an extra approximation for the simulations 



that use Cloudy cooling, but given that most of the free- 
free emission at high temperatures comes from primor- 
dial gas, this is a reasonable approximation. We use n' e , 
instead of the average of the simulated electron number 
density n e , since X-ray observations directly measure the 
value of (n^) from the thermal emission in the cluster, 
not (n e ) (see, e.g., Simionescu et al. 2011). This mea- 
surement therefore produces a bias in the X-ray inferred 
electron density equal to the clumping factor, C, of the 
gas, where 



C 



(6) 



where the averages are volume- weighted. 



3.2.3. 



Central Entropy 

For the purpose of comparing the entropy in the simu- 
lated clusters to the observationally deduced entropy as 
in the H10 sample, we calculate a slightly modified ver- 
sion of the entropy from our simulation grid. While for 
the most part, the variation from the standard thermal 
entropy are small using this method, it is not identical. 
The standard thermal central entropy calculation for the 
simulation thermal properties is 



K 



(kT) 

2/3 ' 



(7) 



where T is the thermal temperature of each simulation 
zone in the central region, and n e is the mass-weighted 
average of the central electron density in those same re- 
gions. When comparing to observational data, however, 
this is not precisely analogous. The deduced tempera- 
ture from X-ray data is the spectroscopic temperature, 
which is not exactly equal to the mass-weighted mean 
gas temperature in our simulations. Therefore, we re- 
place the average of kT in this calculation with the value 
of fcTgpcc, which is the spectroscopic-like temperature of 
Rasia et al. (2005), described in §3.2.1. We also use n' e 
(defined in §3.2.2) in the entropy calculation to account 
for the bias introduced by X-ray determinations of the 
electron density, and allow a clean comparison with ob- 
servations. A more observationally analogous quantity 
to be generated from our simulations is then 



K 



kT K , 



/2/3 



(8) 



over the central region of interest, which we use for our 
comparisons to observations. 

3.2.4. Central Temperature Ratio 

We calculate the central temperature ratio as the ratio 
of the central temperature To = T spoCi o over the virial 
temperature T vir . The virial temperature is calculated 
from a scaling relation of M 50 o (500 x the mean back- 
ground density; Eqn. 4 of H10): 



lkcV 



logi 



(i7(i 



-^500 

2.5 x lO^h^M, 







(9) 



We use this mass scaling relation instead of calculating 
it from the mass or spectroscopic-like weighted temper- 
ature because it eliminates the need to remove cold gas 
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(see §4.5) in some of the clusters (which is equivalent 
to determining a cluster has a cool core before we try 
to categorize it later). It is appropriate for clusters at 
low redshift with masses similar to our simulated clus- 
ters, and it is also how observationally-derived masses 
are calculated in H10. 

4. UNDERSTANDING THE EFFECTS OF VARIED 
NUMERICAL METHODS 

Our first suite of simulations focus on the effects of 
Cloudy cooling (i.e., metal-dependent) and distributed 
feedback on simulations and observables. This is ac- 
complished with four simulations beginning with an ini- 
tial simulation similar to what was performed in B08. 
The only functional difference is a somewhat higher en- 
ergy feedback parameter e_g: 4.11 x 10 -6 in B08 versus 
1.0 x 10 -5 here. Our current value releases an amount 
of energy per unit stellar mass within a factor of two 
of several other similar studies of clusters (e.g. Kravtsov 
ct al. 2005; Tornatorc ct al. 2007). A further difference is 
that in B08 each cluster is simulated in a separate zoom- 
in simulation, while in this work we simulate all of the 
clusters at once. This should not, however, lead to any 
meaningful differences because the mass and peak spatial 
resolutions, refinement criteria, and the physics modules 
are identical in the two cases. 

We label the first simulation "RS-Single" (i.e. metal- 
independent cooling with feedback deposited in a single 
cell) to make clear which cooling and feedback modules 
are in use. The remaining three simulations (see Table 
1) swap RS cooling in favor of Cloudy cooling ("Cloudy- 
Single"), swap single cell for distributed feedback ("RS- 
Dist"), or swap both ( "Cloudy-HZ-LE" , which can be 
read as "Cloudy-High Metal Feedback- Low Energy Feed- 
back"). We do not change anything else between simu- 
lations, and in particular, €e and ez are kept constant 
between simulations. The figures in this section include 
all the gas in the central regions of the clusters, which 
may include cold gas invisible to X-ray observatories. We 
discuss the impact of the cold gas on our calculations in 
§4.5. 

Additionally, below we compare the observable quan- 
tities of our simulated clusters against the observations 
presented in H10. Figure 2 shows that our most mas- 
sive clusters are less massive than the most massive ob- 
served clusters. Therefore, for purposes of comparing 
simulations to observations, a subset of the H10 sample 
has been selected to exclude clusters with M500 values 
greater than the most massive simulated clusters. How- 
ever, as we will discuss in §5.1, more massive simulated 
clusters do not exhibit substantially different behavior 
than the clusters discussed in this section. This selec- 
tion eliminates 3 of 28 SCC clusters, 4 of 18 WCCs, and 
5 of 18 NCCs from the observational sample, but does 
not greatly affect the overall ranges, nor distribution of 
observed central cooling time, entropy, or temperature 
ratio. 

4.1. The Impact of Metallicity-Dependent Cooling 

The left column of panels in Figure 1 visually com- 
pares RS-Single with Cloudy-Single using a representa- 
tive cluster. This particular cluster is used because it 
is undergoing a single merger, which allows us to illus- 
trate several important differences between simulations, 



but it is not overly complicated by multiple simultaneous 
mergers. It is clear that the choice of cooling method has 
a pronounced effect on the cluster. The Cloudy-Single 
cluster is less centrally dense than the same cluster in the 
other three models, and the infalling clump that is very 
cold in the RS-Single case is much warmer with smoother 
gradients in, e.g., Cloudy-Single. These differences are 
also evident in the temperature profiles shown in Fig. 3. 
At all radii, the mean Cloudy-Single cluster is warmer 
than the mean RS-Single cluster. Outside O.I-R500, the 
profiles are shifted by a nearly constant amount (~ 0.1), 
which is a direct result of Cloudy cooling preventing un- 
enriched gas at high redshift from becoming too cold. 
The large differences in temperatures inside O.I-R500 are 
also a result of the early overcooling of the RS method 
due to its assumption of relatively high metallicity at all 
times. 

In Fig. 4 we compare simulations against observations, 
and by extension, to one other. Figure 4 contains a great 
deal of information and therefore warrants a detailed de- 
scription. For a given simulation, there arc three ele- 
ments: the main panel and two histograms. The main 
panel (which has the largest area and contains the sim- 
ulation label) shows the central cooling time versus en- 
tropy for both simulations and observations, with the 
simulated clusters represented by circles and the sub- 
selected H10 cluster sample represented by triangles. 
The H10 triangles are oriented according to the central 
cooling time subgroup cuts (SCC < 1 h -1 / 2 Gyr, 1 < 
WCC < 7.7, and 7.7 < NCC). The colors of the glyphs 
in the main panel correspond to the central temperature 
ratio following the color scale along the top of the figure. 
Both above and to the right of each main panel are nor- 
malized histograms of central cooling time (above) and 
entropy (right) for simulations and H10 data. The colors 
of the histogram bars indicate the type of data plotted 
(simulated or H10 subgroups) following the histogram 
color key, which is above all the figure panels and below 
the color bar. The entropy histograms are cut according 
to SCC < 22 h- 1 / 3 keV cm 2 , 22 < WCC < 150, and 150 
< NCC. Histogram bars are semi-transparent, and over- 
lap between simulation and H10 data are indicated by 
the three darker colors in the histogram color key. Note 
that the entropy distributions for panels (a),(b), and (e) 
are highly peaked and one or more histogram bar extends 
beyond the displayed fractional range (0.25). 

We can now remark on the profound effect Cloudy cool- 
ing has on the observables between panels (a) and (c) of 
Fig. 4. While central cooling times are much shorter 
than observed in the RS-Single case, Cloudy-Single clus- 
ters exhibit a much wider distribution of cooling times 
that are shifted to higher, more physically realistic val- 
ues. A similar effect is seen in the central entropy values 
as well. However, despite the apparent visual improve- 
ment, non-parametric Kolmogorov-Smirnov (KS) (Kol- 
mogorov 1933; Smirnov 1939) tests between data and 
simulation for these two observables gives p-values less 
than 0.01 (the typical lower cutoff for statistical similar- 
ity) for both RS-Single and Cloudy-Single, which indi- 
cates that none of the simulated distributions are statis- 
tically similar to observations. The short cooling times 
and low entropies with RS cooling are a direct result of 
early universal overcooling. The overall colder gas pro- 
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Fig. 1 . — Projections of gas column density (left panel in each pair) and gas density-weighted temperature (right panel) of the corresponding 
halo from the four simulations as labeled. The density and temperature color bars apply to all four pairs. "RS-Single" and "RS-Dist" 
employ metallicity-independent gas cooling, while "Cloudy-Single" and "Cloudy-HZ-LE" use metallicity-dependent cooling. Star particles 
in "RS-Single" and "Cloudy-Single" deposit all their feedback into the single cell that surrounds them, while stars in "RS-Dist" and 
"Cloudy-HZ-LE" use distributed feedback and deposit their feedback over the 27 cells that surround them. The projections are through a 
cube 2 Mpc on a side. The cluster has R200 ~ 1-9 Mpc, M200 ~ 3.9 X 10 14 M Q and T vir PS 3.7 X 10 7 K (using Eqn. 9). 



duces higher central cluster densities (Fig. 5), and lower 
central temperatures (e.g., Fig. 3) which, in the relevant 
temperature range (~ 10 7 K), result in higher cooling 
rates (Eqn. 5). Note that Fig. 5 shows that outside 
the core, both methods of cooling produce cluster den- 
sity profiles in rough agreement with observations, but 
that towards the core all simulations show central densi- 
ties higher than observations, which is a typical problem 
with these kinds of simulations. The use of Cloudy cool- 
ing creates clusters with less unphysically high central 
densities than with RS cooling, but additional physical 
models (e.g. AGN feedback) may be required to address 
this discrepancy fully. 

In Figure 6 we compare the central temperature ra- 
tios of simulated and observed clusters. It's clear that 
neither RS-Single (a) nor Cloudy-Single (c) are similar 
to observations. RS-Single is too cool and produces few 
NCC clusters, while Cloudy-Single is the opposite, with 
too many NCC clusters. Neither of the KS tests be- 
tween these two simulations and the observed distribu- 
tion of central temperature ratios give p-values greater 
than 0.01. 

4.2. The Impact of Distributed Supernova Feedback 

Next, we study the effect of varying the spatial extent 
of thermal and metal feedback from supernovae (using 
RS cooling in both simulations). Figure 7 shows that the 



mean temperature profiles of RS-Single and RS-Dist are 
indistinguishable. This is a reasonable result because the 
addition of distributed feedback does not change the total 
amount of energy deposited into the gas by a given star 
particle. Interestingly, there are some differences that 
can be seen when the temperature projections are com- 
pared in Fig. 1. From RS-Single to RS-Dist, it appears 
that the amount of very cold gas in the infalling clump 
is somewhat diminished. This is consistent with the dis- 
tributed energy feedback of merger-driven star formation 
more effectively heating the cold gas of the clump over a 
wider region. This same effect is seen in projections of 
other merging cold clumps in RS-Dist clusters. Note that 
although although the infalling clump has less cold gas 
with distributed feedback, the mass of the clump is low 
enough (and other merging clumps like it) that averaged 
over all the clusters, the warming effect is not enough 
to dramatically change the mean temperature profile be- 
tween the two simulations. 

The effect of distributed feedback on observables is 
smaller than the changes when Cloudy cooling is used, 
which is illustrated in Fig. 4. The distributions of cen- 
tral cooling times and entropies of the RS-Dist clusters 
are slightly widened, but do not resemble observations. 
Likewise, the spread of temperature ratios (Fig. 6) is 
slightly greater with distributed feedback, but again is 
not very similar to observations. Like RS-Single, none of 
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Fig. 2. — Normalized distribution of M500 for simulated and ob- 
served clusters. The colors are semi-transparent - dark blue indi- 
cates overlap between observation and simulation. The values of 
M500 for H10 data are found using their Eqn. 4 and fcT v ; r from 
their Table 2. The simulations in the lower four panels are intro- 
duced and discussed in Sections 5 and 6. 

the KS tests for RS-Dist indicate that the simulated dis- 
tributions are statistically similar to observations. Con- 
sidering all three observables together, it appears that 
distributed feedback has a slightly similar, but much less 
strong effect on the observables as does Cloudy. Namely, 
by distributing feedback over wider regions, including 
cold, tight clumps (like in Fig. 1) it acts to partially 
counter the RS overcooling and widen the distributions 
of the observables. The entropies for both RS simulations 
are low because nearly all the clusters have cool central 
temperatures and central densities higher than Cloudy 
simulations (Fig. 5). These two effects are both a con- 
sequence of early overcooling, and are complementary in 
Eqn. 8, resulting in low entropies. 
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Fig. 3. — The mean (solid lines) normalized temperature profiles 
for all Cloudy-Single (blue) and RS-Single (red) clusters and the 
la scatter (semi-transparent filled regions). 



4.3. The Impact of Combining Distributed Feedback + 
Metallicity-Dependent Cooling 

Finally, we examine the result when both distributed 
feedback and Cloudy cooling are used simultaneously. 
Morphologically, the combination of Cloudy cooling and 
distributed feedback appears to only make small differ- 
ences when compared to the substitution of only Cloudy. 
Indeed, the qualitative differences between rows of Fig. 
1 are greater than the differences between columns. This 
highlights the importance of the addition of metallicity- 
dependent cooling due to its greater effect on the clus- 
ters when compared to the addition of distributed feed- 
back. This is reasonable because distributed feedback is 
only operational when and where star particles are being 
formed, while Cloudy cooling makes substantial changes 
to cooling rates compared to RS cooling over the whole 
simulation at all times. 

Figure 8 illustrates this point. Outside O.I-R500, which 
accounts for the bulk of the cluster, the temperature 
profiles are nearly identical. However, inside O.I-R500, 
where the star particles are applying feedback, the differ- 
ences between Cloudy-Single and Cloudy-HZ-LE become 
much larger. In particular, the la scatter of Cloudy-HZ- 
LE temperatures is much larger than it is for Cloudy- 
Single. Although Cloudy-Single produces approximately 
25% more star particles by mass than Cloudy-HZ-LE (see 
Fig. 13 in §6.1), the ratio of total mass of metals (in star 
particles and in the gas) to the total mass of stars is the 
same within few percent in both simulations. The dif- 
ference is that in Cloudy-HZ-LE, 19% of the metals by 
mass are stored in the gas, while this figure is only 10% 
for Cloudy- Single, which implies that the Cloudy-HZ-LE 
has ~ 50% more metals by mass in the gas than Cloudy- 
Single. The increased metals in the gas allows for much 
higher rates of cooling, giving rise to the scatter seen in 
Fig. 8. 

The fact that more metals are stored in star particles 
in Cloudy-Single is very important. Although some of 
the metals given to star particles are recycled and en- 
riched back into the gas after the particle is formed (see 
Eqn. 2), the metals that remain in star particles after 
the feedback is applied (e.g., after 12idyn) are locked in 
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Fig. 4. — Comparison of central cooling time, entropy, and temperature ratio for both simulated and observed clusters. Please see text 
for a detailed description of this figure. The circled points in panel (d) correspond to the two clusters in Fig. 9. 

forever. Crucially, these metals are no longer available 
to contribute to cooling of the intracluster medium. The 
net effect is that, in comparison to the Cloudy-HZ-LE 
star particles, the Cloudy-Single star particles behave as 
a kind of "metal sink." 

This result of enhanced metal feedback in the gas is 
a curious secondary effect of distributed feedback. It 
appears that without distributed feedback, the metals 
returned by star particles in each star-forming region stay 
contained within that region, and are re-captured by the 



next generation of star particles. Distributed feedback 
appears to disrupt this cycle because it deposits much 
of the metal-enriched gas outside star formation regions, 
preventing metals from being locked in star particles, and 
resulting in more metal mass in the intracluster medium. 
We compare the metallicities for these two simulations to 
observations in §4.4. 

The impact of the higher metallicity and cooling in 
Cloudy-HZ-LE is starkly visible in Fig. 4(c,d). Instead 
of a tight linear grouping like the Cloudy-Single clusters, 
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Fig. 5. — Mean density profiles for each simulation (solid lines) 
with la scatter (semi-transparent filled regions). The range of 
observed densities from Planck Collaboration et al. (2011) is over- 
plotted in grey. 

the Cloudy-HZ-LE clusters are scattered widely in cen- 
tral cooling time, entropy, and temperature ratio. Many 
clusters exhibit quite low cooling times, entropies, and 
extremely low central temperatures, and some of these 
"cold-core" clusters have T /T vir much less than 0.01 and 
have central temperatures that drop below 10 6 K. The 
cold-cores are tightly bunched at the low end of the dis- 
tribution in Fig. 6(d). It is these cold-cores that give 
rise to the scatter in Cloudy-HZ-LE cluster temperature 
profiles inside 0.1 i?200 (Fig- 8). As with the earlier 
simulations, none of the Cloudy-HZ-LE distributions are 
statistically similar to observations according to the KS 
test. We note that in reality, X-ray observations cannot 
detect the cold gas present in the cold-cores, which we 
discuss in more detail in §4.5. 

Using a halo merger tree we examine the evolution of 
the clusters, and, in particular, focus on the central tem- 
perature histories of the cold-core clusters. As we found 
in B08, after formation all of the clusters in our calcu- 
lations eventually become CC clusters with low central 
temperatures. Later on, some of the clusters experience 
major mergers, which may disrupt the embryonic cool 
core, and raises the central temperatures permanently. 
However, some clusters grow by only minor mergers and 
smooth accretion early on, and the low central temper- 
atures can be preserved to z=0, although they may en- 
counter large impact parameters and/or minor mergers 
at late times. These later mergers are insufficient to de- 
stroy the large and tightly-bound cool core, and only 
act to add more cool gas to the core.' With the addi- 
tion of Cloudy cooling, and the enhanced distribution of 
metals from distributed feedback algorithm, the Cloudy- 
HZ-LE clusters that might have stayed simply cool are 
in fact able to become cold. This is illustrated for two 
Cloudy-HZ-LE clusters in Fig. 9, where the central tem- 
perature ratio is shown along with major merger events 
over time. Note that the mass-scaling relation used to 
calculate T vir elsewhere in this paper is not used in this 
figure because at early times the mass and redshift of 
the clusters falls outside the applicable range of the scal- 
ing relation. Instead, T v j r is calculated by finding the 
spectroscopic- like temperature inside i?200j but with the 
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Central Temperature Ratio (To/T v - ir ) 

Fig. 6. — Central temperature ratio histograms for simulated 
clusters (green) and observed clusters of comparable mass. The 
calculations for the simulated clusters do not include a tempera- 
ture cut. The colors for observed data correspond to the bimodal 
groups NCC (red) and CC (blue) as defined by H10 (CC < 0.7 = 
To /T v i r ). The simulations in the lower four panels are introduced 
and discussed in Sections 5 and 6. 

core removed. At z=0 this T V j r is roughly similar to the 
mass-scaled T v i r . In Fig. 9 the cluster that ends up 
warm at z=0 is cool from zf«3 to zssl.3, but then un- 
dergoes two mergers in quick succession around z«1.3. 
After a period of mixing, over which the core temporar- 
ily drops in temperature, the central temperature rises 
and never re-cools. In contrast, the cold cluster expe- 
riences no major mergers early in its formation, which 
allows the cold core to be established. It does experience 
several mergers later in its evolution, but they happen 
after the cold core has become established and too grav- 
itationally bound to be destroyed by the energy of the 
various mergers. We stress that while this scenario of 
the simulated cool/cold-core clusters surviving late ma- 
jor mergers may be inconsistent with observations that 
show clusters with recent major mergers lack cool cores 
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Fig. 7. — Similar to Fig. 3, except the simulations used are RS- 
Single (blue line) and RS-Dist (green line). 
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Fig. 8. — Similar to Fig. 3, except the simulations used are 
Cloudy-Single (blue line) and Cloudy-HZ-LE (green line). 

(e.g., Rossetti et al. 2011), late mergers impacting CCs 
is not a requirement of our model. 

4.4. Metallicity 

We have just shown that simulations using RS 
(metallicity-independent) cooling fail to reproduce the 
observed central quantities due to unphysically high lev- 
els of cooling early in the cluster evolution. Unfortu- 
nately, using a metal-dependent cooling algorithm does 
not result in realistic clusters either. It is known that a 
bimodality exists between CC and NCC clusters in which 
NCC clusters exhibit a flat metallicity profile towards 
the center while CC clusters show an enhancement in 
the core (De Grandi & Molendi 2001; Bohringer et al. 
2004; Baldi et al. 2007; Johnson et al. 2011). Clearly, 
gas metallicity is a very important quantity to exam- 
ine when attempting to understand CC/NCC clusters, 
and in Fig. 10 we compare the metal profiles of both 
Cloudy simulations to observations (Matsushita 2011). 
These observations, done with the using XMM-Newton 
satellite, estimate metallicity by measuring the strength 
of iron lines as function of radius in 28 clusters, 20 of 
which are also in H10 sample. Fig. 10 mirrors the ear- 
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Fig. 9. — To/T v i r as a function of time for the two Cloudy- 
HZ-LE clusters indicated by circles in Fig. 4(d). The red line 
corresponds to the warm cluster in the upper right, and the blue 
to the cold cluster in the lower left. The downward-pointing arrows 
show major (M su bhalo/M c i ustcr > 0.1) merger events for the warm 
cluster. Likewise, the upward arrows are for the cold cluster. The 
number with each arrow is the M su t,halo/M c i ustcr ratio for that 
merger. At z=0 the warm cluster has M200 ~ 6.0 X 10 14 Mq, and 
the cold cluster has M 2 oo ~ 1-5 x 10 14 M . 



lier finding that the Cloudy-HZ-LE simulation has higher 
mean metallicity than Cloudy-Single, but it also shows 
two important differences between simulations and ob- 
servations. First, the central (< 0.07i?2oo) metallicities 
of both Cloudy simulations are far too high, and second, 
the slopes of the profiles are far too steep with too little 
metallicity outside > 0.07i?2oo- Both effects have been 
noted before in simulations with star formation + cool- 
ing and without AGN feedback (e.g., Sijacki & Springcl 
2006). 
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Fig. 10. — Mean Cloudy-Single and Cloudy-HZ-LE metallicity 
profiles (solid lines) and la scatter (semi-transparent shaded re- 
gions). Observed data from Matsushita (2011) is shown with error 
bars; blue for clusters with a cD galaxy, and orange for those with- 
out. 

Both effects are caused by star formation that is too 
centrally concentrated. Nearly all of the star formation 
occurs inside O.I.R200 m the region of the cluster with 
the highest gas density, which implies that essentially all 
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of the thermal feedback is contained within this region. 
This effect is still present when distributed feedback is 
used, since the region over which feedback takes place 
for each star particle is a cube « 47 hr l kpc on a side 
- significantly smaller than the volume contained within 
0.1i?2oo- Given the high density of the gas, cooling takes 
place rapidly, lowering the entropy of metal-enriched gas 
and keeping it contained within the center of the cluster. 
In principle, star formation in galaxies near the outskirts 
of clusters could result in substantial additional enrich- 
ment of the gas. In practice, however, the relatively low 
resolution of these calculations prohibits star formation 
in all but the largest galaxies in a cluster, resulting in 
highly centralized star formation. As a result of these 
two factors, the metallicity profile in all of our simulated 
galaxy clusters is heavily tilted, with the inner regions 
of the clusters having far too much metal, and the outer 
regions having too little. 

4.5. The Effect of Cold Gas on Observables 

The orbiting X-ray observatories Chandra and XMM- 
Newton have little to no sensitivity below sa 0.5 keV, 
which means that they are unable to accurately detect 
gas colder than a few million K (Garmire et al. 2003; 
Jansen et al. 2001). The simulated cold-core clusters 
discussed earlier in this section contain gas at their cen- 
ter well below this temperature threshold, which means 
that the cold gas would be invisible to the observatories. 
Therefore, to more accurately model observables as they 
might be measured by the aforementioned observatories, 
it is appropriate for the calculations of the observables 
to remove the cold central gas from consideration (Nagai 
& Lau 2011). As described in Section 3.2.1, mock X- 
ray measurements like T spec are unphysical at these low 
temperatures. Calculations of the metallicity and other 
X-ray derived quantities from simulations are likewise 
skewed if this non- X-ray emitting gas is included. 

The removal of cold gas from the calculation is accom- 
plished for each cluster by first finding the total central 
X-ray emissivity in the 0.5-7.0 keV bandpass. Next, we 
find the temperature for which gas below that temper- 
ature accounts for only 5% of the total emissivity. We 
then eliminate the gas below this threshold from our cal- 
culations of all three main observables. The mean tem- 
perature cut across all clusters is roughly 1 — 2 x 10 7 
K, which is very similar to the lower temperature limit 
of the observatories quoted earlier. We show the result 
of these temperature cuts for two simulations, Cloudy- 
Single and Cloudy-HZ-LE, in Fig. 4(e, f). We note that 
when we perform this same temperature cut procedure 
on the other simulations (not plotted), the differences in 
the observables are negligible because those simulations 
contain very little or no gas below 1 — 2 x 10 7 K in the 
centers of the clusters. 

With the temperature cuts, the cold-core clusters dis- 
appear, and the simulations appear more similar to ob- 
servations (see panels (e, f) of Fig. 4). In fact, Cloudy- 
HZ-LE without cold gas produces a distribution of cen- 
tral cooling times just on the cusp of being likely to be 
consistent with observations (p=0.01). However, we ar- 
gue that this improvement is of mixed value for analyzing 
these simulations. While the cold-core clusters do ap- 
pear more realistic using the temperature cuts, it hides 
the true nature of the cold gas which is unlikely to be 



physical. For example, Hicks et al. (2010) find using UV 
GALEX observations of CC clusters (some of which are 
in the H10 sample) that star formation rates in the inner 
50-100 kpc are generally well under 1 M /yr. The low 
rate of star formation in observed CC clusters indicates 
that there is not a large amount of very cold gas present 
in the centers of CC clusters. Because the cold gas is 
therefore unlikely to be physical, and clusters without 
the unphysical cold gas are unaffected by the tempera- 
ture cuts, the inclusion of the cold gas in our calculations 
is actually more informative for purposes of understand- 
ing the simulations and identifying problematic clusters. 
Therefore, we will include the cold gas in our calculations 
for the remainder of our analyses. 

5. THE EFFECT OF VARIED FEEDBACK AND CLUSTER 
MASS ON OBSERVABLES 

The results presented in the previous demonstrate that 
simply changing radiative cooling methods and the way 
that stellar feedback is deposited into the intraclustcr 
medium cannot produce physically-reasonable simulated 
clusters that agree with observations. Therefore, we 
move on to varying the quantities of feedback into the 
ICM from supernovae, which are remaining free param- 
eters in our model. The results from Section 4 indi- 
cate that the RS cooling results (using a metallicity- 
independent cooling table) fail to produce realistic clus- 
ters, which leaves the Cloudy (e.g., metallicity-dependent 
cooling) simulations as options for a platform for fur- 
ther exploration. It is arguable which of the two Cloudy 
cooling simulations better reproduce observations, but 
we choose to use Cloudy-HZ-LE as the basis for further 
testing because it employs both metallicity-dependent 
cooling and the distributed stellar feedback algorithm. 
The Cloudy cooling has no adjustable parameters, leav- 
ing only the feedback parameters e# and ez- 

Because the central metallicity of Cloudy-HZ-LE is so 
much higher than observations, and it is clear that the 
high metallicity has a large effect on the observable prop- 
erties of the simulated cluster cores, we run three simu- 
lations in which we lower the metal feedback parameter 
(tz) by a factor of five, and we vary the value of the en- 
ergy feedback parameter (eg) (see Table 1). These three 
simulations are labeled "Cloudy-LZ-LE" ("Cloudy-Low 
Metal Feedback-Low Energy Feedback"), "Cloudy-LZ- 
ME" (Medium Energy Feedback), and "Cloudy-LZ-HE" 
(High Energy Feedback). The first of the three has the 
same value of as Cloudy-HZ-LE, while the next two 
increase by a factor of 8 and 20, respectively. 

Figure 11 shows that the lower ez has dramatically 
reduced the central metallicities in all three simula- 
tions. Indeed, inside roughly 0.05i?2oo (which includes 
our defined central region) , the Cloudy-LZ-LE metallici- 
ties agree fairly well with observations. However, outside 
that region the metallicity still declines far too rapidly 
for all three simulations. This is not surprising since, to 
first order, we have simply reduced the amount of metal 
that is produced overall in the simulation, but not the 
locations or timing of star formation. 

In Figure 12 we show the results of varying the feed- 
back energy on our observables. The figure shows that 
the lower and more realistic central metallicities have a 
profound effect on the observables. From Cloudy-HZ- 
LE (Fig. 4(d)) to Cloudy-LZ-LE, the very low entropy 
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Fig. 11. — Mean metallicity profiles for all three low metal Cloudy 
runs (solid lines) and lcr scatter (semi-transparent filled regions). 
Observed data from Matsushita (2011) is shown with error bars; 
blue for clusters with a cD galaxy, and orange for those without. 

tail of Cloudy-HZ-LE has vanished, and in its place the 
central entropy distribution is now peaked solidly, but 
too strongly, in the observed SCC range. Additionally, 
lowering the metal feedback has significantly altered the 
distribution of central temperature ratios (Fig. 6(e-g)), 
and brought them closer to agreement with observa- 
tions. In particular, the KS test between observations 
and the Cloudy-LZ-LE central temperature ratio distri- 
bution gives p — 0.03, which is the best match out of all 
the distributions presented in this paper. 

The other two low-metallicity runs show that by ad- 
justing the feedback parameters, we are able to produce 
more realistic clusters, at least as measured by these 
three observables. In particular, except for a few sim- 
ulated cold-core clusters (which we discuss below), the 
agreement between the observed SCC+WCC distribu- 
tions of central cooling time and of entropy and the 
Cloudy-LZ-HE clusters is visually better than any of the 
other simulations. As in all other simulations, KS tests 
produce p- values less than 0.01, but this is largely due to 
our inability to produce clusters with NCC-like cooling 
times and entropies. 

These three low metallicity simulations produce a few 
interesting trends that allow us to highlight the inter- 
connectedness of star formation, gas enrichment, and gas 
cooling that is at the core of our inability to produce a 
realistic distribution of clusters. First, the simple act 
of reducing the metal feedback from Cloudy-HZ-LE to 
Cloudy-LZ-LE noticeably shifts the distributions of cool- 
ing times and entropies to lower values. This is a direct 
result of the fact that lowering the amount of metals in 
the gas reduces gas cooling, and this in turn reduces the 
star formation rate (see Fig. 13). Lower rates of star for- 
mation reduces the amount of energy feedback available 
for raising the central cooling times and entropies. We 
discuss Figure 13 in more detail in Sections 6.1 and 6.2. 

Second, as ce is increased, the bulk of the central cool- 
ing times and entropies moves to greater values, and be- 
comes more in line with the SCC+WCC portion of the 
observed distributions. This is not a surprising effect. 
Raising the energy feedback of star formation acts to 
heat up the gas, which raises cooling time, and also in- 



creases the entropy. Note that despite increasing by 
a factor of 20, we still fail to produce fully-NCC (i.e. in 
the NCC range of all three central quantities) clusters. 

Curiously, increasing the energy feedback also in- 
creases the number of very low entropy, cold-core clus- 
ters from Cloudy-LZ-LE to Cloudy-LZ-HE. In contrast 
to the Cloudy-HZ-LE cold cores, which are primarily due 
to the high metallicity of the ICM in the center of the 
cluster, the increase in the number of cold cores seen 
in the Cloudy-LZ-HE simulation are due to a combina- 
tion of low metallicities and low amounts of formation 
of stars and their associated feedback energy. Figure 14 
shows that the coldest Cloudy-LZ-HE clusters also have 
relatively low total stellar masses. In the higher energy 
feedback run, it appears that the high thermal feedback 
of a given generation of stars inhibits the formation of 
the next generation of stars, forming a feedback loop. 
In some clusters this feedback loop is very strong, and 
despite (and due to) the high value of eg, the star forma- 
tion is not enough to warm up the core. This suppression 
of star formation as feedback energy is increased is also 
illustrated in Fig. 13. 

Third, while the central cooling times and entropies 
shift as the thermal feedback is modified, the central 
temperature ratios (Fig. 6(e-g)) do not change nearly 
as much (except for the addition of the few cold-core 
clusters at the cold end of the distributions). The rea- 
sons why the central temperatures do not change very 
much with eg, while the cooling times and entropies do, 
are complex. Figures 11 and 13 show that the increase 
in tE lowers both the central metallicities and star for- 
mation rates. This results in lower cooling rates which 
prevents central densities from becoming quite as high 
(see 15). Keeping roughly equal temperatures, but hav- 
ing lower densities, results in longer cooling times and 
higher entropies. Across the three simulations the tem- 
peratures stay roughly equal because, despite the fac- 
tor of 20 increase in the amount of thermal energy re- 
turned to the intracluster medium per solar mass of star 
formed, the amount of feedback energy returned to the 
gas stays strikingly constant (Fig 16). From Cloudy-LZ- 
LE to Cloudy-LZ-HE, there is less than a factor of two 
increase in total energy deposited in the gas. Figure 16 
is discussed in more detail in §6.3. 

5.1. Higher Mass Clusters 

Above and in §4 we removed 12 clusters from the H10 
sample because their masses are higher than the most 
massive ones in our simulations. Due to their deeper 
potential wells and likely richer history of mergers, it 
is reasonable to expect that higher-mass clusters would 
tend toward having NCCs when compared to less mas- 
sive clusters. In order to test if higher-mass clusters be- 
have differently than the ones already discussed, we ran 
a suite of follow-up simulations of higher-mass clusters. 
The higher mass clusters do not show any meaningful 
changes in their physical observables when compared to 
their lower mass analogues, indicating that simply pro- 
ducing higher mass clusters does not produce more real- 
istic clusters. 

6. DISCUSSION 

In this paper, we have shown that Cloudy (metallicity- 
dependent) cooling and distributed thermal feedback can 
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Fig. 12. — Plots and histograms of central cooling time and entropy for simulated and observed clusters; similar to Fig. 4. See §4 for 
a detailed description of the aspects of this figure. Gas temperature cuts have not been applied to the calculations in this figure. The 
"Adiabatic" simulation is discussed in Section 6.3. 



be used to create simulated galaxy clusters with some- 
what more realistic ICM properties than in previous gen- 
erations of simulations. Even our best models, however, 
fail to accurately reproduce the observed cluster proper- 
ties - in particular, the distribution of clusters into cool 
core and non-cool core clusters and related properties: 
the central entropy and cooling time of the intracluster 
gas, the ratio in temperature of this central gas from the 
cluster's virial temperature, and the distribution of met- 
als within the cluster itself. In particular, virtually none 
of our simulated clusters fall into the non-cool core ranges 
of central cooling time and entropy. In this section, we 
discuss possible reasons for the continued difference be- 
tween our simulation clusters and real galaxy clusters. 

As discussed in the introduction to this paper, there 
is a large body of numerical work that shows that AGN 
can help to alleviate some of the problems exhibited by 
our simulations, including regulating star formation and 
promoting the distribution of metals. In Section 7 we 
discuss a sub-grid feedback method that we believe will 
allow us to address the some of problems discussed in this 
section and should allow us to include AGN feedback to 
our clusters as well. 



6.1. Lack of Resolution 

Due to finite computational resources and the need to 
simulate a large volume of the Universe, our galaxy clus- 
ter simulations have relatively poor mass and spatial res- 
olution compared to the current cutting-edge in galaxy 
or single-cluster (e.g., Li & Bryan 2012) formation sim- 
ulations. This has a negative effect on the accuracy of 
our model, and is especially true at early times, when 
cooling and star formation is taking place in real galax- 
ies with masses below the resolution of these simulations. 
For example, at z« 6 the dominant galaxy population in 
terms of feedback is w 10 10 — 10 11 M (Bouwens et al. 
2007; Gonzalez et al. 2011), and at z=2 it is « 10 12 M Q 
(Conroy & Wechsler 2009). Given that the dark matter 
particle mass in our simulations is 7.8 x 10 9 M Q , we sim- 
ply cannot resolve galaxies at high z, and therefore we do 
not make stars at that epoch. This is clearly illustrated 
in Fig. 13. Furthermore, in Enzo, as in hydrodynami- 
cal simulations in general, gas cooling within halos only 
becomes significant once once the gravitational potential 
is well-resolved, and thus can compress gas to high den- 
sities and form stars. From a practical standpoint, this 
means that a halo must contain a relatively large number 
of dark matter particles before star formation can occur 
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Fig. 13. — Star formation rates over time for each simulation. 
We plot observed rates as a visual guide, but we do not expect 
concordance. The simulated rates have been smoothed by averag- 
ing over a moving window ~ 200 Myr wide. The blue left-pointing 
triangles are observed rates from Bouwens et al. (2007), and the 
green right-pointing are from Hopkins (2004, 2007). Error bars 
have been omitted from the observations for clarity. 
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Fig. 15. — Mean density profiles for each simulation as labeled 
(solid lines) with la scatter (semi-transparent filled regions). The 
range of observed densities from Planck Collaboration ct al. (2011) 
is over-plotted in grey. These density profiles are essentially un- 
changed from the earlier simulations using Cloudy cooling (Fig. 
5). 
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- substantially more than is required to simply register 
the existence of a halo. This is illustrated by the fact 
that in these simulations the average halo has a mass 
of 7.3 x 10 12 M© when a star is first formed inside it, 
which corresponds towl, 000 dark matter particles. Be- 
cause the halos that form stars for the first time are so 
massive, star formation occurs later in time and in ob- 
jects with deeper potential wells (and thus denser gas) 
than is typical in the real Universe. The denser gas is 
able to cool more readily, which means that the thermal 
feedback from star formation does not warm up gas as 
effectively, and the overall result is an excess of cooling 
and a deficit of entropy at low redshifts. The injection of 
metals into the intergalactic medium is similarly affected 
by our lack of resolution, with metal-enriched gas being 
more strongly concentrated than is observed. This leads 
to an unfortunate feedback loop - stars form too late and 
in over-dense regions, which are then over-polluted with 
metals. This excess of metal in the gas causes cooling to 
occur too rapidly, reducing or eliminating the pressure 
gradients that would drive metal-enriched gas to large 
cluster radii. 

The star formation rate in our simulations is illustrated 
in Figure 13. Observed data from Bouwcns et al. (2007) 
(from the HST HUDF and GOODS fields) and Hopkins 
(2004, 2007) (combined from the literature) are included 
in the figure for comparison. This illustrates that all of 
the calculations discussed in this paper differ substan- 
tially from observed rates at z>0.3, and in some cases 
do so quite dramatically. Of course, as mentioned above, 
the simulations do not resolve the galaxies that are re- 
sponsible for a large fraction of star formation, and there- 
fore it is unreasonable to expect precisely matching rates. 
Rather, the figure confirms that star formation begins 
later than the cosmic mean, which is exactly the oppo- 
site of what one would expect from observations of galaxy 
clusters, where the bulk of stars form earlier than the 
cosmic average due to clusters being substantially over- 
dense regions at high redshift. 

6.2. Interconnectedness of Star Formation and Feedback 

A striking feature of the different stellar histories in 
Fig. 13 is how much they change as the models for ra- 
diative cooling and subgrid star formation are adjusted. 
Every deliberate choice we made to improve the distribu- 
tions of our three core observables has, in fact, caused the 
timing of the peak rate of star formation to occur later 
and later. Recall that this is contrary to what is desired, 
which is a cluster star formation rate peak that occurs 
earlier than the overall universal peak. This effect, and 
the effects noted earlier of adjusting ce and ez, illustrate 
the interconnectedness of star formation and feedback in 
our simulations. 

An example of this effect are the low-metallicity runs 
discussed in §5. The first run, Cloudy-LZ-LE, that only 
reduces the metal feedback, results in clusters that have 
systematically low central entropies and cooling times. 
Increasing the energy feedback does somewhat amelio- 
rate those particular discrepancies, but the result is an 
even less realistic history of star formation. It is clear 
that the star formation — > metal feedback —¥ gas cooling 
cycle is a self-reinforcing cycle. A logical change to one 
aspect of our subgrid stellar model results in other as- 
pects drifting further from the observational results we 



are attempting to match. As mentioned previously, the 
cold-core clusters seen in Cloudy-LZ-HE are evidence of 
an incorrectly modeled cycle of star formation and metal 
enrichment. This same behavior exists in all the other 
simulations discussed in this work, and it is impossible 
to disentangle a change in one part of that cycle without 
affecting all the others, often negatively. 

More broadly, the fundamental challenge with our sub- 
grid models for star formation and feedback is their 
strong dependence on the mass and spatial resolution 
of the simulation, as was discussed in Section 6.1. Using 
our current technique, we will not be able to correctly 
resolve the early evolution of the galaxy populations in 
proto-cluster regions, which limits our ability to bring 
the observable properties of our simulated clusters in line 
with real galaxy clusters. Given this unresolvable issue, 
it is clear that a fundamentally different approach needs 
to be taken, which we will elaborate upon below. 

6.3. Energetics 

Our lack of fully-NCC clusters indicate that we are not 
correctly applying thermal feedback to these simulations. 
As an experiment, we run a simulation that has identical 
initial conditions and simulation parameters to the other 
calculations described in this paper, except that radia- 
tive cooling, star formation, and stellar feedback are all 
disabled. We label this calculation "Adiabatic," and the 
results are shown in Figures 6(h) and 12(d). The only 
source of cooling in this calculation is from adiabatic ex- 
pansion of the plasma, and the only sources of heating 
comes from adiabatic compression and from shock heat- 
ing, powered by gravitational potential energy. As can be 
seen in Figures 6(h) and 12(d), the Adiabatic simulation 
results in warm galaxy cluster cores, with the majority 
of simulated clusters solidly in the NCC range of the cen- 
tral temperature ratio distribution (Fig. 6(h)). However, 
the central cooling times (as before, calculated using a 
cooling curve based on an optically thin plasma with a 
metallicity of 0.5) and entropies are only in the WCC 
ranges (Fig. 12(d)), which shows that simply heating the 
core to NCC-like temperatures is not enough to produce 
fully-NCC clusters. The Adiabatic results also highlight 
the requirement that clusters include both heating and 
cooling because without the early cool phase of clusters, 
some of which become may NCC clusters following merg- 
ers, it is impossible to create the full spectrum of cluster 
types. 

We now turn to Figure 16, which demonstrates the in- 
sensitivity of the observables in our simulations to our 
choice of cooling algorithm and magnitude of both metal 
and thermal feedback. The instantaneous energies are 
calculated for each simulation by convolving the star for- 
mation rate (Fig. 13) with Eqn. 1 and the value of ee 
(Table 1). Following the rates in the top panel, the maxi- 
mum instantaneous ratio between any two simulations is 
under an order of magnitude. The bottom panel is per- 
haps even more interesting. This shows that across all 
seven simulations the total feedback energies deposited 
in the gas are remarkably similar. From the highest to- 
tal in RS-Single, to the lowest in Cloudy-LZ-LE, there 
is only about a factor of three change in the amount of 
energy deposited, which stands in stark contrast to the 
large changes in physics modules and feedback parame- 
ters. This shows that using our current tools it is very 
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difficult to change the energy feedback history, and that 
it is likely very challenging to simulate the required en- 
ergy feedback history to produce realistic clusters with 
the parameters we are able to adjust. 

Despite the assertions above about little energetic dif- 
ferences between simulations, there are significant dif- 
ferences in the observables across the simulations. The 
fact that even small changes in energetics lead to major 
shifts in (especially) entropy and cooling time highlights 
the importance of accurately modeling energy physics in 
cluster simulations. As mentioned before, we discuss in 
§7 a method that should allow us to more accurately 
model the energetics in clusters and reproduce more re- 
alistic results. 

7. CONCLUSIONS AND FUTURE WORK 

Where are we along the road to making realistic galaxy 
clusters? In Table 2 we summarize qualitatively the suc- 
cesses and failures of our seven main simulations, basing 
our statements on the data contained in Figures 4, 6, 
and 12. Since no simulations produced fully-non-cool 
clusters (as defined by the measured cooling times and 
entropies), we omit that as a comparison from the table, 
and note that it is a failure of all the simulations. How- 
ever, we also note that all simulations, except those using 
metal-independent RS cooling, did produce a substantial 
fraction of NCC clusters as measured by central temper- 
ature ratio. Where we feel the judgment is warranted, 
we have labelled assessments with either a "best" or a 
"worst." Our assessments of the three main observables 
illustrate that overall none of the simulations did par- 
ticularly well at reproducing the observed distributions. 
The simulations using cooling tables that assume a con- 
stant, high metallicity (the RS series of runs) are the 
furthest from matching observations, primarily due to 
over-cooling of gas. Among runs using metal-dependent 
cooling (the Cloudy series of runs), the simulation with 
the low metal feedback and high thermal energy feed- 
back (Cloudy-LZ-HE) overall produces the best match to 
observations, with reasonable agreement for both cool- 
ing time and decent agreement for both cooling times 
and entropies, and central temperatures almost entirely 
within the observed range, although the full spectrum of 
cool-core and non-cool core cluster is not reproduced. 

We have also included assessments of the metallicities 
and star formation histories of the simulations in Table 
2. In short, they generally fail to agree with observa- 
tions. In terms of metallicity, Cloudy-LZ-LE is the best 
because it has central metallicities that are in relatively 
good agreement with observations. That said, this cal- 
culation, like all the other simulations, shows too steep 
a gradient in the metal profile at larger radii. With re- 
gards to star formation rate, none of our simulations ac- 
curately capture the star formation history one expects 
from galaxy clusters: the overall rate of star formation 
is too low, and peaks at too late of a time, compared to 
the rates inferred from observations of cluster galaxies. 

Notwithstanding the qualitative assessments above 
and in the table, some of our results are encouraging. 
The inclusion of metallicity-dependent cooling rates sig- 
nificantly improves upon previous work that uses a cool- 
ing table that assumes a fixed, relatively high metallic- 
ity, and we show that it is possible to adjust simula- 
tions to achieve better specific agreement with observa- 



tions. Clearly, a metal-dependent cooling method, like 
the method of Smith et al. (2008, 2011), is essential to 
prevent early overcooling of unenriched gas. Using a dis- 
tributed feedback model, rather than depositing all of the 
thermal energy from supernovae into a single cell, does 
not have as large of an effect on cluster properties as cool- 
ing. In addition, we find that while metallicity profiles as 
a function of radius can be altered somewhat based on 
choices of feedback parameters, the metallicity profiles 
never resemble cluster observations. This points to the 
need for a more sophisticated feedback model that may 
include AGN feedback. Taken as a whole, these results 
begin to indicate a way forward that uses metallicity- 
dependent cooling, and which also addresses the short- 
comings of the star formation and feedback model em- 
ployed in this paper. 

We conclude that our current methodologies are not 
capable of reproducing all observed properties of real 
clusters. Therefore, the next step on the road to more 
realistic clusters is to employ new physical models that 
address the shortcomings discussed in this paper. In par- 
ticular, it is clear that we must address the unrealistic 
star formation histories and metallicity profiles of the 
clusters. Without improvements to both, it will be im- 
possible to simulate the energetics of heating and cooling 
of observed clusters. 

The method we believe will provide substantial im- 
provements to our current methods are "galaxy par- 
ticles" (also called "galaxy constructs", or "galcons," 
Arieli et al. 2008, 2010). Galaxy particles employ a model 
where the the time evolution of mass, metal, and en- 
ergy feedback of galaxies are derived from the observed 
globally-averaged star formation rate, and the feedback 
quantities are deposited into the intergalactic medium 
isotropically and over a sphere whose radius is substan- 
tially larger than the cell size of the simulation, thus 
potentially avoiding the over-cooling. Arieli et al. (2010) 
compares a cluster simulation employing galaxy parti- 
cles to one similar to RS-Single (in particular, in terms 
of cooling and feedback models). They find that galaxy 
particles are more effective at depositing mass, energy, 
and metals into the intracluster medium which results 
in more realistic temperature, density, metallicity, and 
entropy profiles. They find that galaxy particles help re- 
duce the amount of subhalo "over-merging" in clusters, 
which reduces the central dominance of the simulated 
cluster. An additional advantage of such a technique 
is that it will allow us to include, in a parametric way, 
more complex physics such as AGN feedback and merger- 
driven starburst events. If successful, one could then use 
a wider range of cluster observables (e.g., the Sunyaev- 
Zel'dovich effect or L^-T relationship) to further con- 
strain calculations. Galaxy particles are presently being 
worked on by the Enzo development community, and we 
will present the results of cluster simulations using galaxy 
particles in a future paper. 
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TABLE 2 

Simulation Successes and Failures 



Simulation Label Cooling Time 



Entropy 



Temperature Ratio 



Mctallicity Profile 



Star Formation 



Mostly CCs, 
Too Few NCCs 

Mostly CCs, 
Too Few NCCs 
Too Many NCCs, 
Too Few CCs 

Cloudy-HZ-LE Decent SCC+WCC, Many Cold-Cores With Too Bimodal 



RS-Singlc Worst, Entirely 
Too Short 
RS-Dist Mainly Too Short 

Cloudy-Single Best SCC+WCC 
agreement 



Mainly Too Low 

Mainly Too Low 

Too Peaked Between 
SCCs and WCCs 



But A Few Too Low 



Too Low Values 



N/A 

N/A 

High in Center, 
Profile Too Steep 
Worst, Far Too High in 
Center, Profile Too Steep 



Best Overall, Too Low 
Early, Too Slow Rolloff 
Too Low Early, 
Too Slow Rolloff 
Too Low Early, 
Too Slow Rolloff 
Far Too Low Early, 
Too Slow Rolloff 



Cloudy-LZ-LE Generally Too Short Mainly SCCs, Peaked Between the Best in Center, Low in Far Too Low Early, 

A Few WCCs CC & NCC Groups Periphery, Profile Too Steep Too Slow Rolloff 

Cloudy-LZ-ME Mainly SCC, Almost Entirely Peaked Between the Low in Periphery, Far Too Low Early, 

A Few WCCs WCCs CC & NCC Groups Profile Too Steep No Rolloff 

Cloudy-LZ-HE Decent SCC+WCC, Decent SCC+WCC, Peaked Between the Low Overall, Worst, Abysmal Always, 

A Few Too Short A Few Too Low CC and NCC Groups Profile Too Steep No Rolloff 

Note. — The entries are qualitative assessments of how well the simulations agree with observations. 
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